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Marginal quantization of an Euler diffusion process and its 

application 



Abass Sagna 



Abstract 



We propose a new approach to quantize the marginals of the discrete Euler process resulting from the 
discretization of a brownian diffusion process using the Euler scheme. The method is built recursively 
using the distribution of the marginals of the discrete Euler process. The quantization error associated to 
the marginals is shown to goes toward at the optimal rate associated to the quantization of an M'^-valued 
random vector In the one dimensional setting we illustrate how to perform the optimal grids using the 
0^ Newton algorithm and show how to estimate the associated weights from a recursive formula. Numerical 

^^ tests are carried out for the pricing of European options in a local volatility model and a comparison 

^> with the Monte Carlo simulations shows that the proposed method is more efficient than the Monte Carlo 

n , method. 

cd 1 Introduction 

a 

'~~' Optimal quantization method appears first in EOll where the author studies in particular the optimal quanti- 

^-H zation problem for the uniform distribution. It has become an important field of information theory since the 

K^ early 1940's. A common use of quantization is the conversion of a continuous signal into a discrete signal 

_^ that assumes only a finite number of values. 

l/^ Since then, optimal quantization is applied in many fields as in Physics, in Computer Sciences, . . . , 

^ and recently in Numerical Probability from the seminal work |[T2l . Its application to Numerical Probability 

■^ relies on the possibility to resume a random vector taking values in a set of infinite cardinality by a discrete 

^ random vector valued in a set of finite cardinality with the associated weights. This allows to approximate 

,-H either expectations or more significantly some conditional expectations. This is the key property used to 

^jj solve some problems arising in Quantitative Finance as optimal stopping problems (see |[I]|21), the pricing 

of swing options (see Q), stochastic control problems (see EO), nonlinear filtering problems (see e.g. 
llT3l[T7l l5l[T8l). the pricing of barrier options (see [ 19|). 
C^ In Quantitative Finance, several problems of interest are the estimation of quantities like (for a given 

function f -.W^^R) 

E[/(Xt)], T>0, (1) 

or involving terms like 

E[f{Xt)\Xs = x], 0<s<t, (2) 

where {Xt)te[o,T] is a stochastic process which evolves following the stochastic differential equation 

dXt = bit, Xt)dt + a{t, Xt)dWu ^0 = so G M"', (3) 

where VF is a standard d-dimensional Brownian motion starting at and where the functions 6 : [0, T] x 
^d |_^ ^d ^jjj jj^g matrix diffusion coefficient function a : [0, T] x M'^ i-> M^ x M"' are measurable and satisfy 



^Federation de Mathematiques d'Evry, Laboratoire Analyse et Probabilites, 23 Boulevard de France, 91037 Evry, & ENSUE, 
e-mail: abass. sagna@ensiie.fr. This research benefited from the support of the " Chaire Marches en Mutation", Federation 
Bancaire Frangaise. 



some appropriate conditions which ensure the existence of a strong solution of the stochastic differential 
equation. Since in general the solution of ([3]) is not explicit we have first to approximate the continuous 
paths of the process (-^i)igro,T] by discrete paths using some discretization scheme like the Euler scheme. 
Given the (regular) time discretization steps tk = A;A, A; = 0, . . . , n, A = T /n, the "discrete Euler process" 
{Xt^)k, k = 0, . . . ,n, associated to the previous diffusion process (Xi)^gro,T] i^ defined recursively as 

Xt,^, = Xt,+b{tk,Xt,)A + a{tk,Xt,){Wt,^, - Wt,), Xo = Xq. 

Then, once we have access to the discrete paths of the stochastic process {Xt)t£[o,T]' the quantities (|T| and 
(|2]) are estimated by 

mixtj] (4) 

and 

E[f{Xt,+,)\Xt, = x], tk+i =t,tk = s. (5) 

Remark 1.1. a) Under some assumptions on / (for e.g. if / is four times continuously differentiable with 
polynomial growth) the estimation of E(/(Xr) by E(/(Xt)) induces the following weak error convergence: 

\Ef{XT)-Ef{XT)\<- 

n 

with C > and where n is the number of time discretization steps, 
b) Also remark that for every p>l, 

e( sup \Xtj\ <+oo. (6) 

The estimation of quantities like Q or (|5]l can be performed using Monte Carlo simulations. Neverthe- 
less, an alternative to the Monte Carlo method may be the optimal quantization method, specially in small 
dimension (d < 4 in the theory but in practice it may stay competitive with respect to the Monte Carlo 
method up to dimension d = 10, see |[T6l ). 

In fact, suppose that we have access to the optimal quantization or to some "good" (in a sense to be 
specified later) quantizations {Xf^ '' )fe of the process {Xt^.)k on the grids x^* = {x-^ *,•••, Xj^/'} of size 

Nf., for k = 0, . . . ,n. Suppose also that we can get the associated weights ^{X^^ '^ = x^ *), z = 1, . . . , A^^, 

k = 0, . . . ,n, and the transition probabilities p-'/^ix) = ^{X'f^ '^ = x,- ''*^\Xf^ ^ = x) for every k = 
0, . . . , n — 1. Then using optimal quantization method, the expressions Q and ([5]) are estimated respectively 
by 



and 






The question is then to know how to get the optimal grids x^*=, for A: = 0, . . . , n, the associated weights 
and transition probabilities. In a more general framework, as soon as the stochastic process {Xt^)k (or the 
underlying diffusion process {Xt)t>o) can be simulated one may use stochastic (or Lloyd) algorithms to 
estimate the (optimal) grids and the associated weights or transition probabilities. In the special case of the 
one dimensional setting we can use the Newton algorithm in some situations. This deterministic algorithm 
leads to more precise estimations and is more fast (at far in some cases) than stochastic algorithms. 

Consider the first quantity Q, given that the second one is of interest in a working paper. Suppose for 
example that we aim to estimate the price of a Put option with a maturity T, a strike K, an interest rate r and 



a given present value Xq in tiie Black-Sciioles model wiiere tlie dynamics of the stock price process under 
the risk neutral probability is given by: 

dXt = rXtdt + aXtdWt, Xq = xq G M, 

which solution reads 

Xt = Xo e(^— V2)t+^m 4 Xo e(--V2)t+av^z^ ^ _ _^(0. ^y 

Then we want to estimate the quantity (keep in mind that this is a toy example given that it can be written in 
a semi-closed formula involving the cumulative distribution function of the standard Gaussian distribution: 
the famous Black-Scholes formula) 

e-'^E{f{XT)), 

where the payoff function /(x) = m.ax{K — x, 0). To estimate such an expression by optimal quantization, 
one may quantize the random variable Z and compute the associated weights using the Newton algorithm. 
Then, setting Z^ its optimal quantization on the grid x^ = {x^ , . . . , x^} we estimate the quantity of 
interest by 

i=l 

where g : x ^ Xq e^'''~" /2)T+crvTx ^^^^ ^ g [j2]). Now, suppose that the model is a local volatility model 
where the dynamics of the stock price process evolves following the stochastic differential equation (called 
Pseudo-CEV in QOJ): 

dXt = rXtdt + {} ^ ' dWt, Xo = xo (7) 

for some 6 G (0, 1) and ?? G (0, ??] with t9 > 0, where r is the interest rate. In this situation the distribution 
of Xt is not known and if we want to estimate the quantity of interest: e^'^^E(/(Xr), where f{x) := 
max(/f — x,0) is the payoff function, we have to discretize the process {Xt)t>o as {Xtf,)k=o,...,n, with 
tn = T, using for e.g. the Euler scheme, then, estimate 

e-'-^E(/(XT)) 

by optimal quantization. Up to now, the only way to get the optimal grids and the associated weights in this 
situation is carrying out stochastic algorithms (see e.g. JTH) or Lloyd's algorithms (see e.g. [8|), even in the 
one dimensional framework. As pointed out previously this method may be very time consuming. 

In this paper we propose a new approach to quantize the process {Xtf.)k, based on a recursive method 
involving the conditional distribution of the marginals: means, those of Xt,,^^\Xtf,. This approach raises 
some interesting problems among with the computation of the rate of convergence of the quantization error 
associated to Xt^ (for /c = 0, . . . , n), the question of knowing if the defined quantization Xf^ of Xt^ on 
the grid x^'' satisfy the stationary property, . . . 

We observe that by construction these marginal quantizations are stationary. This property is very useful 
when estimating expressions like ([4]) and (|5]l by optimal quantization. We also observe that the rate of 
convergence of the quantization error associated to the marginal Xt,, is of O^N^"^), where N]^ is the level 
of the quantization. As a consequence of this result, if / is a Lipschitz function then estimating E/(Xt) by 
Ef{X^") where X^ " is the optimal quantization of Xt of size Nn will induce the following error bound: 

|E/(XT)-E/(X^")| = 0(iV-i). 

On the other hand, if A; G {0, . . . , n — 1}, we show that one may rely the L^-quantization error E|Xt^^^ — 
Xf """^ P associated to Xt, . , with the distribution of Xt, as 

nX,^^, - Xll'^' |2 = /" e( . min \Y,ix) - xf +^ hp(Xt, G dx) (8) 



where Yk{x) is the d-dimensional gaussian distribution with known mean and variance. Thus, in practice, to 
quantize the random variable Xt^._^_^ we have to compute an expectation with respect to Xt^ ■ The distribution 
of Xt^. is not expUcit but since the random variable Xt^ i^ supposed to be already quantized we will estimate 
(HJby 

5^e( . min |n(xf'=) - xf +^nP(xf '= = xf'^) (9) 

where Xf^ ^ is the quantization of Xt^. on the grid x^^ . We show that the error induced by such an estimation 

iVr. 

is bounded by the quantization error ||Xjj. — Xf ||2. 

From the numerical viewpoint our method allows us to estimate efficiently, as well from the precision 
of the estimates as from the computational time, expressions like Q. In fact, in the specific case of the 
one dimensional setting, this approach allows us to use the Newton algorithm as an alternative to stochastic 
algorithms. This is a very important fact because it is more precise and more fast than stochastic algorithms. 
Notice however that when dealing with high dimension {d > 2) there is at the moment no way to avoid the 
use of stochastic (or Lloyd) like algorithms. 

The paper is organized as follows. We recall first some basic facts about optimal quantization. The 
marginal quantization method is described in Section [3] We give in this section the induced quantization 
error and compute the error bound derived from the estimation of ([8]) by Q. We also illustrate how to get the 
optimal grids using Newton's algorithm and show how to estimate the associated weights. The last section. 
Section [4j is devoted to the application of the proposed method to the pricing of an European Put option in 
a local volatility model (as well as in the Black-Scholes model) and to the comparison with the Monte Carlo 
method. 

2 Overview on optimal quantization methods 

Let {fl, A, P) be a probability space and let X : {Q, A, P) — > M'^ be a random variable with distribution 
Fx- The L''(Px) -optimal quantization problem of size N for the random vector X (or for the distribution 
Fx) consists in finding the best approximation of X by a Borel function of X taking at most N values. 
Assuming that X G L^{F), we associated to every Borel function q{X) taking at most A^ values the U- 
error ||X — g(X)||r measuring the distance between the two random vectors X and q{X) w.r.t. the U- 
norm, where \\X\\^ := (E|X|^) '^ and where | • | denotes an arbitrary norm on R"^. Then finding the best 
approximation of X by a Borel function of X taking at most A^ values turns out to find the solution of the 
following minimization problem: 

cnA^) =ini{\\X -q{X)\\r,q:R'^^x^,x^ C M'^, card(x^) < N}. 

Now, let x^ = {x^ , . . . , x'^} C M'^ be a subset (a codebook) of size N (it is also called an A'^-quantizer 
or a grid of size N) and define a Voronoi partition Ci(x^)j=i at of M*^, which is a Borel partition of M*^ 
satisfying for every i € {!,..., N}, 



Ci(x^)c{xeR'^:\x-xf\= min |x-xf|}. 

7=1,. ..,A' ^ 



Consider the Voronoi quantization of X (or simply called quantization of X) on the A^-quantizer x^ defined 
by 

TV 
vx'^ _ V^ N-. 

i=l 

Then for any Borel function g : M*^ — )• x^ = {x'^ , . . . , x^} we have 



\X-q{X)\> mm d{X,xf^) = d{X,x^) = \X - X'''' \ P a.s 



so that the optimal quantization error eN,r{X) reads 



-2:^11 ^N ^ ju,d „„^At„N\ 



eN,r{^) = inf{||X-X^ ||^,x''' CM'', card(x''') <iV} 

^ l/r 
1 / ay^,X^ 

card(a;^)<Ar 



inf / d{z,x'^Y(mx{z) . (10) 



Note that for every A^ > 1, the infimum in (10) is reached at one codebook at least. Any A^-quantizer 
realizing this infimum is called an ^''-optimal A^-quantizer Morever, if card(supp(Px)) > N then the 
optimal A^-quantizer is of size exactly N {see ^ or ifTlll ). On the other hand, the quantization error, eN,r{X), 
decreases to zero as the grid size A^ goes to infinity and its rate of convergence is ruled by the so-called Zador 
Theorem given below. 

Theorem 2.1. fZador Theorem, see f9']) : Let X be an W^ -valued random vector such that ¥.\X\^~^^ < 
+00 for some ij > and let ¥x = -Pa + -Ps b^ the Lebesgue decomposition of¥x with respect to the 
Lebesgue measure Xd, where Pa denotes the absolutely continuous part and Ps the singular part. Then 

lim N''/''{eNAP)y = Qri^x) (11) 

with 

d+r 

Qr{^x) = Jr,d( I f^dXd] ' =Jr,d 11/11^ G[0,+Oo), 

Jr4 = mf^ A^'-/'^e^,,([/([0, 1]'^)) e (0, +oo), 

where C/([0, 1] ) denotes the uniform distribution on the set [0, 1] and f = -^. 

We will call Qri^x) the Zador's constant associated to X. From the Numerical Probability point of 
view, finding an optimal A^-quantizer x^ may be a challenging task. In practice (we will only consider the 
quadratic case, i.e. when r = 2 for numerical implementations) we are sometimes led to find some "good" 
quantizations X^ which are close to X in distribution, so that for every Borel function / : M i— )■ M, we 
can approximate E[/(X)] by 

N 

E/(l-")=^/(xf)p„ (12) 

where pi = P(X^ = xf). Among "good" quantizations of X we have stationary quantizers defined as 
follows. 

Definition 2.1. An A^-quantizer x^ = {x^, . . . , x^} inducing the quantization X^ of X is said stationary 
if 

V i / j, xi / Xj and P (AT G Ui5C,(x^)) = 

and 

e[a:|x^'^J =1^'^. (13) 

Defining the distortion function by 

DM,2ix'')= [ d{z,x''fdFx{z) = y2 [ \z-xf\^dFx{z), (14) 

a stationary quantizer x^ = {x^ , . . . ,x^} is in fact an A^-quantizer satisfying the stationary equality: 
^Dn,2{x^) = 0. The following result justify the introversion of the differentiation and the integral leading 
to ([T3]l when differentiating ([14|, see m. 



Proposition 2.2. The function Dj\^^2 i^ dijferentiable at any N-tuple x G (M ) having pairwise distinct 
components and satisfying 

¥{X eUidC^{x^)) =0 

(which condition holds in particular when X is a continuous random vector) and we have 

VZ?Ar,2(x^) = (2 / [x^ -z)dFx{z)], . (15) 



C,{xN) 



=1,...,N 



For numerical implementations, the stationary quantizers search is based on zero search recursive pro- 
cedures like Newton's algorithm in the one dimensional framework and some algorithms as Lloyd's I al- 
gorithms (see e.g. |8]), the Competitive Learning Vector Quantization (CLVQ) algorithm (see El) or 
stochastic algorithms (see iflSl ) in the multidimensional framework. One may download optimal quan- 
tizers associated to some random variables as the multivariate Gaussian random vector in the website 
www . quantize .math-f i . com. 

When estimating the expectation Ef{X) by E/(X^ ) where x^ is an A^-quantizer and X^ the quan- 
tization of X on x^, we make an error which may be bounded by the quantization error K\X — X^ p, 
depending to the regularity of the function /. We next recall some error bounds induced from the approxi- 



mation of Kf{X) by ( 12 1, we refer to 1 16|) for further details. 

a) Let x^ be a stationary quantizer and / be a Borel function on R*^. If / is a convex function then 

E/(X-'')<E/(X). (16) 

b) Lipschitz functions: 

- If / is Lipschitz continuous then for any A^-quantizer x^ we have 

|E/(X)-E/(X-")| < [/]Lip||X-X-"||2, 

where 

\fix)-fiy)\ 
[/J Lip := sup , _ , . 

Xy^y 1-^ y\ 

- Let : M'^ — > M+ be a nonnegative convex function such that 9{X) G L^(P). If / is locally 
Lipschitz with at most ^-growth, i.e. \f{x) - f{y)\ < [/]Lipk - y\{0{x) + 9{y)) then f{X) G 
Li(P) and 

|E/(X) - E/(X-")| < 2[/]Lip||X - X-''U9iX)h. 

c) Differentiable functionals: if / is differentiable on M*^ with an a-Holder differential D/ (a G [0, 1]), 
then for any stationary A^-quantizer x^, 

\EfiX)-Ef{X^'')\ < [Df]4X-X'-''\\l+''. 

3 Marginal quantization of a discrete stochastic diffusion process 

Let {Xt)t>o be a stochastic process taking values in d-dimensional Euclidean space W^ and evolving follow- 
ing the stochastic differential equation: 

dXt = bit, Xt)dt + a{t, Xt)dWt, Xo = xo G M"', (17) 



where 14^ is a standard d-dimensional Brownian motion starting at and where the functions b : [0,T] x 
M*^ I— )■ W^ and the matrix diffusion coefficient function a : [0, T] x M"' i— ;• M'^ x M^ are measurable and 
satisfy the global Lipschitz and linear growth conditions: 

\b{t,x)-b{t,y)\ + \a{t,x)-a{t,y)\ < K\x - y\ (18) 

\b{t,x)\ + \a{t,x)\ < K{l + \x\). (19) 

This guaranties the existence of a strong solution of ([TTJI. We also suppose that the matrix a is positive 
definite. In the rest of the paper we will suppose that M'^ is equipped with the Euclidean norm. 

Let tk = kA, k = 0, 1, . . . , n be time discretization steps, with A > 0. Consider the Euler scheme of 
the process {Xt)t>o starting from Xq = xq: 

Xt,+. = K + Ktk. Xt,)A + a{tk,XtJ(Wt,^^ - Wt,). 

Define the functions ruk and Sfc by 

mi.{x) = X + Ab{tk,x) and T,j^{x) = Aa{tk,x)a{tk,x) . 

It is clear that 

Xtk+i\Xtk = X '^J\f{mk{x);T.k{x)) 

and will denote by ^mk{x),T.k{x){') its density and by Tm^{x),'Ek{x)i') its cumulative distribution function: 
where |Sfc(x)| denote the determinant of Tik{x) and for Xk+i G M'^, 

'^mk{x),T,kix)i^k+l) = / ^mk{x),T,k{x)(y)dy- 

J]~00,Xk + l] 

For every n > 0, we will denote by Fx{Xt^ G dxn) the density of Xt^ when starting the process at 

Xq = X. 

3.1 Distortion function of the marginals 

To compute the distortion function of the marginals of the Euler diffusion process we need to know the 
probability distribution of the marginals. We recall how to compute the distribution of the marginal Xt^,_^^ 
with respect to the distribution of Xt^ . 
Using Fubini's theorem we have 

F{Xt,^^ < xu+i) = I I nXt,^, G dy\Xt^ = Xfc)P(X,, G dxu) 

J]-oo,Xk+i\ JIRd 
= / ^v^k{x^),J^^{xk){xk+l)HXtu^dXk). 



Then 



(20) 



dXk+l ^rrH:(xt,),^t,(x^){Xk+l)F{Xtk ^ dXk) (21) 

k 

/ ••• / Tl'^mi{x,),s,{x,){xi+i)fi{dxo)dxi...dxk (22) 



dxk+1 

t Wtirl f im^ 

i=0 

7 



It follows from Equation ( 2 1 1 that the probability density function of the Euler diffusion process may be 
computed recursively if the distribution of Xq (supposing that Xq is random) is known. 

On the other hand this result allows us to compute the (quadratic) quantization error of the marginal 
distribution ^tj._^j given the distribution of Xj^. and then to deduce from a recursive procedure the optimal 
quantization of Xt^^-^ given the distribution of Xt^. (or, from the numerical point of view, given the optimal 
quantization of Xt^. and its associated weights). In fact, let X§ be the quadratic optimal quantization of 
Xq over the optimal grid x^" = {x^ °, ■ • • , Xj^"} (if Xq = x is not random we set A'^o = 1 and Xq ° = 



x). We will deduce the quantization of the marginal random variable Xt^ using Equation (21 ). Once this 
quantization is perform we will get from the same way the quantization of Xt2 given the distribution of Xt-^ , 
and so on. In the general setting suppose that the quadratic optimal quantizations Xf ,l = 0,...,k,of the 
Xti s over the optimal A'^; -quantizers x'^' = {x^ ',■■■, x^' } are given. We propose in the following, a way 
to deduce the quadratic optimal quantization of Xt^^^^ . 

Let A;G{0,...,n — 1} and let Dfc+i {x^''+^ ) be the distortion function associated to the A'^fe+i-quantizer 



k+1 = {xf '=+^ , . . . , Xj^l+l). Then it follows from (21 ) that 



= / (E/ . [xf^'-zf^^,i.),,:,(.){z)dz)nK^dx), 



SO that we have the following result. 
Proposition 3.1. We have 



Dk+i{x^''+')= [ e( mill \Ykix)-xf''+'f)F{Xtk(^dx), (23) 

where Ife(x) is the d- dimensional gaussian distribution with mean mk{x) and variance Sfc(x). 

Keep in mind that the L*" -quantization error for the marginal random variable Xt,,^-^ is defined as 

eN,^,AK+^) = inf {||^t,+i - Xf^^^" t, x^'=+^ C M^ card(x^'=+i) < iV,+i}. 



A natural question arising after our constructive approach is to compute the rate of convergence of the 
quantization error associated to the marginals. 

3.2 Error analysis 

We observe that for every time step t^, the quantization error goes to at the optimal rate A^^ where N^ 



is the grid size of the optimal A'^fc-quantizer x^* = {x^ *",•••, x^*"}. In fact, it follows from Equation (22 1 
that Xt^ admits a density function for every k =, 1, . . . , n. Furthermore, we deduce from (pll that for every 
r > 1 and for every k = 1, . . . , n, E|X(j.|'' < +oo. Then Zador's theorem applies and we have for every 
A: = 1, . . . , n, for any L^'-optimal quantizer x^* for Xt^, 

lim Nl"e^^^AXt,) = Qr{F^) (24) 

where e^ r("^tfc) •" ^(l^^tfc ~ Xf^ *'Y') is the quantization error associated to the L''-optimal quantizer 
x^* and Qr (P^^ ) is the Zador constant associated to the distribution of Xt^. . 



As mentioning before, the computation of the quadratic optimal quantizers for Xt,.^-^ suppose that we 
have access to the distribution of X^^. . In practice, since we suppose that Xt,, has already be quantized and 
that we have access to its associated weights, we will estimate the distortion function -Dfc+i associated to 
Xtf.^^ and given in (23 1 by 






^fc ^ ^Nk^ 



(25) 



We give in the following an error bound induced by such an approximation. 
Theorem 3.2, Suppose that there is is a non null polynomial function 7^a(a, h) of degree a G N of the form 

(26) 



such that for every x,z € 



ra{a,b) = Y^ Cija'V, Qj >0, i,j en 

i,j<Q 



sup \aa {t,x) — aa {t, z)\ < \x — z\Va{\x\,\z\ 
ie[o,T] 



(27) 



Let n> 1. Then, for any sequence of stationary Nk-quantizer x *" for Xt^, s.t. for every A; = 0, , 



,n. 



\Xt, - XI '11^0, as Nk -^ +00 



we have 

|Dfc+i(x^'=+i) - i)fc+i(x^'=+i)| < K \\Xt, - Xf^' II2, fork = l,...,n-l. 
for some real constant K > {). 



(28) 



Before dealing with the proof of the theorem let us make some remarks on Assumption ( 26 1. 

Remark 3.1. Note that Assumption (26 1 is satisfied when aa'^ is uniformly Lipschitz in space, uniformly 
in time: there is L > s.t. 

|cro-^(t,x) -acr^(t,2;)| 
sup 1 1 < L. 



We also note that if a is symmetric then Assumption (26 1 holds with Pi (a, h) = K'^{2 + a + h), where K is 
the constant appearing in ([TSJ) and ([19]). In fact, if a is symmetric then 



|cro-^(t,x) -crcr^(t,z)| = \a'^ {t , x) - a"^ {t , z)\ 

< \a{t,x)-a{t,z)\{\a{t,x)\ + \a{t,z)\) 

< K^\x- z\{2 + \x\ + \z\), 

where the first inequality follows from the Holder and Minkowski inequalities and the last inequality from 



( 18 1 and ( 19 1. As a consequence. Assumption (26 1 is satisfied by a large class of model in the one dimen- 
sional framework, including the Black-Scholes model and the Pseudo-CEV model. 

Let us prove the theorem now. 



Proof. We have 

4+1 



i?t+i(x^^+0-^fc+i(x^'=+^; 



Ef{Xt,)-Ef{Xl ')\<E\f{Xt,)-f{Xf^ 



x"k-~ 



where 

fix) = e( min \Yk{x) - xf^+^|') = e( . min fj{Yk{x))), 

with fj {Yk{x)) = \Yk{x) - 2;^^'=+' 1^ It follows that for every z G R"^, 

/(x)<e(. min |/,(n(x))-/,(n(z))|)+E(. min /.(^(z))). 

or, in other words 

/(x)-/(z)<e(. min \f,{Ykix))-f,{Ykiz))\). 
\j=l,...,Nk+i / 

Using a symmetric reasoning gives, for every x,z G R"', 

|/(x)-/(z)|<e(. min |/,(n(x))-/,(n(z))|). 

Since for every a,b & M'^, it holds the inequality |ap — |6p < |a — 6|(|a| + |6|), we have for every 

j = l,...,iVfc+i, 

IfjiYkix)) - MYkiz))\ < \Ykix) - Yk{z)\{\Ykix)\ + \Yk{z)\ + 2|xf +^|). 
Using Holder's inequality and the inequality (|a| + |6|)^ < 2(|ap + |5p) yields, for every x, z G M*^ 

\f{x)-f{z)\ < (E\Yk{x)-Yk{z)\'y^'(E\Yk{x)f+E\Yk{z)\^ + A min |xf'=+^|2^ ^/^ 



Now, using the uniform Lipschitz assumption on b{t, x) and Assumption (27 1 leads to 

E\Yk{x)-Ykiz)\^ = E\mkix)-mkiz) + {i:k{x)-^kiz))zf 

< Ki (\x - z\^ + |6(t, x) - bit, z)p + E| (Sfc(x) - Sfe(z))Z|2) 



< Ki|x-Z|2p2a(k|,|z|) 

for some positive and generic real constant Ki which may vary from line to line and where V2a is a poly- 



nomial function of degree 2a of the form p6|). Moreover, using the linear growth assumption ( |T9| ) on the 
coefficient b of the diffusion and assumption ( 27 1 we get, for every x G M'^, 

E\Ykix)\^ = E\mkix) + ^kix)Z\'' < 2(|mfe(x)|2 + E|Sfc(x)Z|2) 

< 2(|mfc(x)|2 + |Sfc(x)|2E|Z|2) 

< K{\mk{x)\^ + \a{tk,x)f) 

< iv:(i + |x| + |x|^) 

where K is a generic positive constant depending in particular on A. 

On the other hand, since E(d(Xij.^^, x^^+i)^) — )• as Nk+i — )• +00, the following asymptotic density 
property of (x^^+i ) in the support of the distribution of Xt^.^-^ holds : 

Ve>0, VxGsupp(P^^ ), 3n,,^,fcGN, VArfc+i>n,,^,fc, 5(x,e) n x^^+i / 0. (29) 

'■fc + 1 



N, 



(p) 



Otherwise, there exists x G supp(P^ ), e > and a subsequence (x '=+i)p>i so that Vp > 1, i?(x, e) n 



A^. 



(p) 



"^fc+i 



X fc+i = 0. Then, for every p > 1, 
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which contradicts the fact that ||(i(Xt^^-^, a^^'°+^)||2 — ^ as N^+i —5- +cxd. 
Then, it follows from (|29]) that there is C > such that 

mm \x- + < C. 

j=l,...,Nk+i ^ 

Combining the previous results gives for every x, z G M"^ 

\fix) - fiz)\ < K\X - Z\({1 + |X| + \Z\ + |X|2 + |z|2)P2a(k|, \z\)) = K\x - z\Vp{\xl \z\f'^, 



for a positive real constant K and for Vjs, a polynomial function of degree /3 of the form ( 26 1. So, using once 
again Holder inequality yields 

Now, we know that for every r G N, E|Xtj.|'' < +oo. Moreover, since Xf is stationary, it follows from 



the previous statement and from inequality ( 16 1 that for every r G N, 



Hence 

Finally we have shown that 



K\f{Xt,)-fixC)\ < ^P*.-^r 



for some positive real constant K. D 

We focus now on the numerical computation of the quadratic optimal quantizers of the marginal random 
variable Xt^^-^ given the probability distribution function of Xt^. Such a task requires the use of some 
algorithms like stochastic algorithms, Lloyd's algorithms (both requiring the computation of the gradient of 
the distortion function) or Newton's algorithm (specially for the one-dimensional setting) which involves the 
gradient and the hessian matrix of the distortion (we refer to lfT6l for more details). 



3.3 How to get optimal quantizers of the marginals? 



Owing to Proposition 2.2 and Equation (22i), the distortion Dfc_|_i(x^*+i) is continuously differentiable as 
a function of the A'^/j+i-quantizer x^''+^ (having pairwise distinct components) and its gradient is given by 



VDfe+i(x^'=+i) = hi (xf +^ - y)P(X,,^, G rfy) 



2 . 



ixf'+' - y)F{Ykix) G dy)F{Xt, G dx)) 

k + l) -^ / j = l,...,Nk + l 



where Yk{x) is the d-dimensional gaussian distribution with mean mk{x) and variance Tik{x). Keep in mind 
that the quadratic optimal quantizer at level Nf^^i (that we still denote by x^^+O of Xt^_^^ is defined as 



^Nk+i g argmin{Dfc+i(z), z C W, card(z) < iVfc+i}. (30) 
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Remark 3.2. a) If x^^+i is a quadratic optimal A'^fc4_i-quantizer of Xt^^-^ and if Xf^ '_^^ denotes the quan- 
tization of Xt^^^ over the grid x^'=+i. Then x^''+'^ is a stationary quantizer: it is such that 

VA.+i(x^'^+i) = 0. 
Equivalently, we have for every j = 1, . . . , N^^i, 



J^,P{Yt,ix) G C7,(x^'=+i))P(Xi, G dx) 



<'^' = ~^\ \:.: . vrrr:: ..^l^ : , . (31) 



(X,,^, GC,-(x^'=+i)) 



where as previously yfc(x) is the d-dimensional gaussian distribution with mean mk{x) and variance Tik{x). 
This also means that 

e{x,,^,\xCt) = xC: 



"fc+i 



b) In the d-dimensional (d > 2) setting, Equation (31 1 allows us to compute stationary quantizers for ^tj._,_i 
given the distribution of Xj^, using Lloyd's type algorithms. 

As mentioned in |[T6]| . the usual stochastic or Lloyd's companion algorithms become quickly intractable 
when the dimension d of the random vector Xt^ is greater or equal to 2 due to the fact that we have to 
compute d-dimension integrals on Voronoi cells. Moreover, in our setting the complexity of the algorithms 
will increase in hight dimension since we have to compute additional d-dimensional integrals. For these 
reasons, we will restrict our analysis to the one-dimensional setting where we will use the Newton algorithm 
to perform recursively quadratic optimal quantizers of the marginals Xt^ given the distribution of Xq. 

3.3.1 The Newton algorithm 

Let VDfc+i (x^'^+i ) and V^Dk-{-i{x'^'=+'^) denote respectively the gradient vector and the hessian matrix 
of the distortion function. Using the Newton algorithm, a zero of the gradient may be computed via the 
following recursive procedure starting from a given j;^fc+i'0 g M^'=+i : 

So we have to compute the gradient and the inverse of the hessian of the distortion function. To simplify 
notations set 

iVfc+i iV,+i Nk+i Nk+i 

■^j-1/2 ~ 2 '■^j+l/2~ 2 ' ^^^'^•^1/2 ~ °°^-^Nk+i+l/2~ ^°°' 

and let 

x^^+l-mk{x) x^^+l^-nikix) 
Xk+i,j-{x) := and Xk+ij+{x) := , 

where Vk{x) = Aa'^(tk, x). The components of the gradient can be computed after some elementary com- 
putations and we have for every j = 1, . . . , N^^i, 

dDk+i{x^''+^) 



dxj'+' 



{xj '=+' - nikix)) \To^i{xk+i,j+ix)) - ro,i(xfc+ij_(x))j 

+ Vk{x){^o,iixk+ij+{x)) - $o,i(a;fc+ij_(x)))}p(Xt, G dx). 
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The diagonal terms of the hessian matrix are given by: 



Q2^fk+i 



^o,i{xk+i,j+{x)) - To,i{xk+ij_{x)) 



1 



^^^a>o,i(x,^i,+(x))(<-^-xf-) 



Avkix) 



^OA{xk+lJ-ix)){x 



Nk+1 _ Nk+-_ 
Xj_i 



¥{Xt^ e dx). 



The sub-diagonal terms are 



m 



Nt. 



dx'p+'dx'.^t' 



R Vk 



^(xf +^ - xf^r)^o,i{^k+i,j-i^)nXt, e dx) 



and the super-diagonals are 

d^Dk+iix^''+^) 1 



dxf'+'dxj!^^' 



Vk{x) 



[X 






$o,i(a;fc+ij+(x))P(Xi, G dx). 



Note that for practical implementations we have to estimate the integral terms with respect to the dis- 
tribution of Xtk appearing in the gradient and the hessian of the distortion function. Since our procedure 
is recursive and we have supposed that Xt^. has already been quantized and that the associated weights: 
IP(Xtj. G Ci{x^^)), i = 1, . . . , A'^fc are accessible, we may estimate them by optimal quantization. The 
Newton procedure (|32]) will then be modified as 



^Affc+l,n+l ^ ^. 



where the components of the modified gradient Dk+i {x^'^+^ ) are given for every j = 1, . . . , Nk+i by 



5Dfc+i(x^'=+i) 



dx:'~ 



=f 



^Y.{i^?^' - mkix^n) (ro,i(x,+i,,+(£f'=)) - ToA£k+i,-ix^n)) 

i=l 



+ Vk{x^')i<^o,i{ik+i,M^^')) - $o,i(xfe+i,,-(^^)))}lP(^i. e aix^')). 
The diagonal terms of the hessian matrix V^I)a,-i-i(2;^'=+^'") are given by: 



Q2^Nk+^ 



Yl [ro,l(^A:+l,j+(xf'=)) - ro,l(xfc+ij_(xf'=)) 



4 = 1 



Avkixf"") 



^ci>o,i(xfe+i,,.(i^))(<'=^^ -xf^-r)lp(x,, G c.(x^^)) 



4ffc(Xi '') 



and its sub-diagonal terms are 

d^bk+ijx^^^^) ^ _i '^'^■ 

The super-diagonals terms are 



E^^(<'^' -<-r)'^o-i(^^+i.:'--(^^))ip(^*^ ^ ^^(^'^'=))- 



a^^fc+i(x^^+i) ^ _i ^'= 

axf-axf-^ ~^irivk{x^^)'^'^' 



Nk , 



AT^^+i ^Nt,+ 



x^ -^)$0,l(:£fc+l,,+ (£f'=))P(^t. G C,(x^^)). 
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The i^'^+i's (and its companions x,^^^, x,^^^, ...) are the modified A^^+i -quantizers induced by the 

approximation of the distribution of the Xtf. 's by the distribution of the Xf^ 's in the computation of the 
gradient and the hessian matrix of the Newton algorithm steps. 

A similar idea combining (vector or functional) optimal quantization with Newton-Raphson zero search 
procedure is used in Q in variance reduction context as an alternative and robust method to simulation 
based recursive importance sampling procedure to estimate the optimal change of measure. Furthermore, 
the convergence of the modified Newton algorithm to the optimal quantizer is shown in the framework of 
171 to be bounded by the quantization error. However, the tools used to show it do not apply directly to our 
context and the investigation of the convergence of our modified Newton algorithm to an optimal quantizer 
is a open question. To perform it, we have to show, among others, that the optimal quantizers of marginals 
are unique. This seems to be not obvious to show even if we may reasonably think that the unicity holds. 

Once we have access to the quadratic optimal quantizers x^*- of the marginals Xt^, for A; = 0, . . . , n 
(which are estimated using the Newton algorithm described previously) we have to compute the associated 
weights ^{Xt^ G Cj{x^'')), j = 1 . . . , Nk, for k = 0, . . . ,n. We show in the next result how to estimate 
them. 

Proposition 3.3. Let x^'^+i be the quadratic optimal quantizer of the marginal random variable Xt^^^. 
Given the quadratic optimal quantizer x^'= of Xt,. and the associated weights F(^t^ G Ci{x^'')), i = 



1, . . . , Nk, the probability P(X(j.^^ G Cj{x^''+^)) is approximated for every j = 1, . . . , N^+i by 



F{xt,^, G c,(x^'^-+^)) ^Y.{^oAxk+ij+{x^n) - 

ro,i(xfc+i,,_(xf'=)))p(Xt, G Ciix^")). (34) 

Proof. For every fc G {1, . . . , n — 1} and for every j = 1, . . . , N^^i, we have 

Nk 

F{Xtk+^ G c,(x^^+i)) = ^p(Xt,^, G c,(x^^+0;^t. e aix^')) (35) 

and 

P(X,,^, G C,(x^^+i);X,, G Q(x^'=))=P(X,, < xf^lJ^;Xt^ G Q(x^^)) 

-P (Xi,+i < xf_\+;2 ; Xt, G a (x^'= )) . (36) 

Now we have for every z G M, 

P(Xi,^, < z; Xt, G a(x^^)) =r ( I P(^t,+i G dx\Xt, = y)F{Xt^ G dy) 

Jc'iix'^k) 



which gives the announced result replacing the quantity F(^Xt^^-^ G Cj{x^''+^); Xt^ G Ci{x'^'')) in (35 1 by 



its approximation in ( 36 1. D 



In the next section we propose an application of the proposed method to the pricing of European options 
in a local volatility model. 
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4 Pricing of European options in a local volatility model 

4.1 The model 

We consider a pseudo-CEV model (see e.g. lITOl ) where the dynamics of the stock price process is ruled by 
the following SDE (under the risk neutral probability) 



dXt = rXtdt + ^- 



X, 



5+1 



dWt, Xq = xq 



for some 5 G (0, 1) and -d S (0,i?] with ?? > 0. The parameter r stands for the interest rate and a{x 



(37) 



^ 



\^l+x 



corresponds to the local volatility function. This model becomes very close to the CEV model, 



specially when the initial value of the stock process Xq is large enough. In this case the local volatility 

o-(x) f« 'dx^'"'. 





i9 


MC (10''') 


MC (10") 


MQ 




0.5 


0.0022 


0.0018 


0.0017 


CI 




[0.0017; 0.0028] 


[0.0017; 0.0019] 






0.6 


0.0113 


0.0111 


0.0110 


CI 




[0.0101; 0.0125] 


[0.0107; 0.0115] 






0.7 


0.0377 


0.0373 


0.0370 


CI 




[0.0353; 0.0401] 


[0.0366; 0.0381] 






0.8 


0.0883 


0.0876 


0.0871 


CI 




[0.0843; 0.0923] 


[0.0863; 0.0886] 






0.9 


0.1696 


0.1659 


0.1649 


CI 




[0.1635; 0.1756] 


[0.1640; 0.1678] 






1.0 


0.267 


0.271 


0.271 


CI 




[0.259; 0.275] 


[0.269; 0.274] 






2.0 


2.423 


2.433 


2.426 


CI 




[2.387; 2.459] 


[2.422; 2.445] 






3.0 


5.424 


5.492 


5.478 


CI 




[5.424; 5.512] 


[5.471; 5.512] 






4.0 


8.893 


8.806 


8.808 


CI 




[8.801; 8.986] 


[8.777; 8.835] 





Table 1 : (Pseudo-CEV model) Comparison of the Put prices obtained from Monte Carlo (MC) simulations (followed 
by the size of the MC in brackets) with associated confidence intervals (CI) and from the marginal quantization (MQ) 
method. The parameters are: r = 0.15; 8 = 0.5; n = 120; iVfe = 400, V/c = 1, . . . , n; T = 1; /iT = 100; Xq = 100; 
and for varying values of -Q. 



We aim to compute the price of a European Put option with payoff {K — Xy)+ = max(/f — Xy, 0), 
where K corresponds to the strike of the option and T to its maturity. Then we have to approximate the 
quantity 

e-'^^E(i^ - Xt)+ 

where E stands for the expectation under the risk neutral probability. If the process {Xtf,)k denotes the 
discrete Euler process at regular time discretization steps tk, with = to < • • • < ^n = ^. associated to the 
diffusion process {Xt)t>o, this turns out to estimate 

e-''^E{K - Xt)+ 

by optimal quantization. We estimate this quantity by the marginal quantization method proposed in this 
paper and compare the numerical results to those obtained from Monte Carlo simulations. 
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K 


MC (10'') 


MC (10'') 


MC (lOO 


MQ 




100 


08.89 


08.81 


08.81 


08.81 


CI 




[08.80; 08.99] 


[08.78; 08.84] 


[08.80; 08.82] 






105 


10.61 


10.60 


10.59 


10.59 


CI 




[10.51; 10.72] 


[10.57; 10.63] 


[10.58; 10.60] 






110 


12.53 


12.57 


12.57 


12.57 


CI 




[12.42; 12.64] 


[12.53; 12.60] 


[12.56; 12.58] 






115 


14.72 


14.74 


14.75 


14.75 


CI 




[14.60; 14.84] 


[14.70; 14.78] 


[14.75; 14.77] 






120 


17.18 


17.10 


17.13 


17.12 


CI 




[17.04; 17.31] 


[17.06; 17.15] 


[17.11; 17.14] 






125 


19.64 


19.69 


19.67 


19.67 


CI 




[19.50; 19.78] 


[19.64; 19.73] 


[19.65; 19.68] 






130 


22.41 


22.32 


22.40 


22.40 


CI 




[22.26; 22.56] 


[22.32; 22.41] 


[22.38; 22.41] 





Table 2: (Pseudo-CEV model) Comparison of the Put prices obtained from Monte Carlo (MC) simulations (followed 
by the size of the MC in brackets) with associated confidence intervals (CI) and from the marginal quantization (MQ) 
method. The pai-ameters are: r = 0.15; n = 120; Nk = 400, Vfc = 1, . . . , n; T = 1; i? = 4; Xq = 100; and for 
varying values of K. 



4.1.1 Numerical results 

To deal with numerical examples we set 6 = 0.5, Xq = 100, and choose the interest rate r = 0.15. We 
discretize the price process using the Euler scheme with n = 120 (regular) discretization steps and quantize 
the Euler marginal processes by our proposed method. We put all the marginal quantization grid sizes N^ 
equals to 400 except for Xq" = Xq = 100 which grid size is Nq = 1. We estimate the price of the Put 
option by 

Nn 

E[{K - XC^] = Y^iK - xf")+P(xf" = xf") (38) 

i=l 

where f„ = T, and where x^" = {x-^ ",•••, Xj^"} is the optimal quantizer of size A^^^ computed from the 
Newton algorithm (with 5 iterations) and where the associated weight are estimated from (|34]). 

We compare the prices obtained from the marginal quantization (MQ) method with those obtained by 
the Monte Carlo (MC) simulations even for various values of -d with a fixed strike K = 100 (see Table [T]) 
or for varying the values of the strike K with a fixed t9 = 4 (see Table |2]l. For the Monte Carlo simulations 
we set the sample size M equal to 10^ and 10^ for K = 100 and we set it equal to 10^, 10^ and 10^ when 
making varying the strike K. 

Remark 4.1. (On the computation time) a) Remark that all the quantization grids x^*" of sizes N}^ = 400, 
for every k = 1, . . . ,n = 120, and there associated weights are obtained in about 1 minute from the Newton 
algorithm with 5 iterations. Computations are performed using Scilab software on a CPU 2.7 GHz and 4 Go 
memory computer. 

b) With the same grid sizes as in a), all the optimal grids and the associated weights may be obtained quite 
instantaneously from the Newton algorithm using the C (or C++) programming language. In this case we 
can compute directly the inverse of the hessian matrix using the results in ||2T1|. 

c) It is clear that once the grids and the associated weights are available the estimation of the price by MQ 
method using the sum ( [38] ) is instantanuous. 

Remark 4.2. (Initialization of the Newton algorithm) Let = to < • • • < ^n be the time discretization steps, 
let Xq = X be the present value of the stock price process and suppose that the grid sizes Nk are equal. Since 
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the random variable Xti ~ J^{rno{x);vQ{x)), in order to compute the (optimal) A'^i -quantizer for Xt-^ we 
initialize the algorithm to vo{x)z^^ + mo{x), where z^^ is the optimal A^^i -quantizer of the M{0; 1). Once 
we get the optimal A'^i -quantization x^^ for Xt^ and its associated weights, we initialize the algorithm to 
x^i to perform the optimal A'^2-quantizer for Xt^ and its associated weights, . . . , and so on, until we get 
the optimal A^„-quantizer for Xt^ and the associated weights. Notice that doing so we observe no failure of 
convergence in all the considered examples. 

Remark 4.3. We show in Figure [T] and Figure |2] two graphics where we depict in the abscissa axis the 
optimal grids (of sizes Nf^ = 150) and in ordinate axis the corresponding weights. The dynamics of the 
price process in Figure [T] is given by 

dXt = rXtdt + aXtdWt, Xq = 86.3 

with r = 0.03, a = 0.05 whereas its dynamics in Figure|2]is given by 

dXt = rXtdt + ^ , * dWt, Xq = 100 

^/i+xi 

with r = 0.15, ?? = 0.7, 5 = 0.5. We observe that the forms of the distributions are close to the Gaussian 
distribution as indicated by the Gaussian concentration bounds for some processes (see e.g. lITTil ). This 
suggests that the Gaussian concentration results may exist for this model, at least for some choice of the 
parameter 6. 

Let us come back to the estimation of the Put price. Since the payoff function f{x) = {K — x)+ is 
Lipschitz and since the quantization X^^ " is stationary, it follows from Zador Theorem that 

E/(X,J-E/(Xf;")|=0(iV-i). 
Keeping in mind that the estimation error induced by the MC method satisfies 

M 

M 

i=l 

for a sample X^ , . . . , X^'^ of size M, we deduce that the MQ method provides a precision as good as the 
MC method when Nn ~ \/M- For our numerical examples, we remark first that in all examples the prices 
obtained by MQ stay in the confidence interval induced by the MC price estimates. On the other hand the 
prices obtained by the MQ method are more precise (more specifically when ?? = 4 and K increasing away 
from 100) than those obtained by the MC method when the sample size M equals 10^ or 10^. Consequently, 
the MQ method seems to be more efficient than the MC when the sample size is less than 10^. However, 
when increasing the sample size to M = 10^ the two prices because closer (up to 10~^); despite everything 
the confidence interval associated to the MC method stays large. 

Remark 4.4. We remark that when the Monte Carlo sample sample size M = 10^ it takes about 2 minutes 
and 30 seconds to get a price using the C programming language on the same computer described previously. 
Then, in this situation, it takes more time to obtain a price by MC method than carrying it out by MQ. 
Furthermore, the MQ method would become much more fast when programming the Newton algorithm 
using the C programming language instead of the Scilab software. 

To strengthen the previous conclusions related to the local volatility model we compare the two methods 
in the Black-Scholes framework where the stock price process evolves following the dynamics: 

dXt = rXtdt + aXtdWt, Xq = 100. 

In this setting the true prices are available and will serve us as the support of comparison. The parameters are 
chosen such that the model is close to the Pseudo-CEV model: r = 0.15 and a ss t9Xq~^. Numerical results 
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\^f(XtJ-^f2f^^tJ =0{M-y^) 



are printed in Tables [3] and Table |4] and heighten our conclusions on the Pseudo-CEV model. We notice that 
in the Black-Scholes method, the estimated prices from the MQ method are close to the true prices (the best 
absolute error is of order 10^^ for a volatility a = 5% and the worse absolute error 2.10"^ is achieved in 
high volatility: a = 40%). This show the robustness of the MQ method even for reasonably high values of 
the volatility. 





O" 


MC (10") 


MCCIO'') 


MQ 


Tme price 


Abs. error 




0.05 


0.0015 


0.00178 


0.00176 


0.00177 


1.10"^ 


CI 




[0.0012; 0.0019] 


[0.0017; 0.0019] 










0.06 


0.0116 


0.0109 


0.0109 


0.0112 


3.10^ 


CI 




[0.0104; 0.0128] 


[0.0106; 0.0113] 










0.07 


0.0365 


0.0370 


0.0369 


0.0373 


4.10-4 


CI 




[0.0342; 0.0387] 


[0.0363; 0.0378] 










0.08 


0.0876 


0.0876 


0.0869 


0.0875 


6.10 * 


CI 




[0.0836; 0.0915] 


[0.0863; 0.0888] 










0.09 


0.1666 


0.1644 


0.1647 


0.1654 


7.10-* 


CI 




[0.1607; 0.1724] 


[0.1622; 0.1658] 










0.10 


0.269 


0.271 


0.271 


0.272 


1.10-^ 


CI 




[0.261; 0.277] 


[0.271; 0.273] 










0.20 


2.444 


2.431 


2.424 


2.427 


3.10-^ 


CI 




[2.410; 2.479] 


[2.420; 2.442] 










0.30 


5.455 


5.469 


5.470 


5.474 


4.10-^ 


CI 




[5.395; 5.515] 


[5.450; 5.549] 










0.40 


8.680 


8.787 


8.790 


8.792 


2.10-^ 


CI 




[8.598; 8.763] 


[8.760; 8.813] 









Table 3: (Black-Scholes model) Comparison of the Put prices obtained from Monte Carlo (MC) simulations (fol- 
lowed by the size of the MC in brackets) with associated confidence intervals (CI) and from the marginal quantization 
(MQ) method with the associated absolute error (Abs. error) w. r. t. the true price. The parameters are: r — 0.15; 
n = 120; Nk = 400, V/c = 1, . . . , n; T = 1; iiT = 100; Xq = 100; and for varying values of a. 
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K 


MC (10^) 


MC (10") 


MQ 


Trae price 


Abs. error 




100 


8.680 


8.787 


8.790 


8.792 


2.10^ 


CI 




[8.598; 8.763] 


[8.760; 8.813] 










105 


10.805 


10.739 


10.744 


10.750 


6.10^ 


CI 




[10.71; 10.90] 


[10.71; 10.90] 










110 


12.86 


12.89 


12.90 


12.91 


1.10^ 


CI 




[12.76; 12.96] 


[12.86; 12.93] 










115 


15.29 


15.24 


15.26 


15.27 


1.10^ 


CI 




[15.18; 15.40] 


[15.21; 15.28] 










120 


17.66 


17.81 


17.79 


17.81 


1.10^ 


CI 




[17.54; 17.79] 


[17.78; 17.85] 










125 


20.56 


20.50 


20.50 


20.52 


1.10^ 


CI 




[20.43; 20.69] 


[20.46; 20.54] 










130 


23.28 


23.37 


23.37 


23.39 


2.10^ 


CI 




[23.14; 23.42] 


[23.34; 23.43] 









Table 4: (Black-Scholes model) Comparison of the Put prices obtained from Monte Carlo (MC) simulations (fol- 
lowed by the size of the MC in brackets) with associated confidence intervals (CI) and from the marginal quantization 
(MQ) method with the associated absolute error (Abs. error) w. r. t. the true price. The parameters are: r — 0.15; 
n = 120; Nk = 400, V/s = 1, . . . , n; T = 1; cr = 40%; Xq = 100; and for vai-ying values of K. 



X,,„ = Xi, + b{X,^A + s{Xi,}^Z. Z - m 1); ^^) = rx: ff(i) = dj 




A',. = 4, lt=kA,k = l 2.'): A = 0.03; ? = l Nt 



Figure 1: ("Black-Scholes model ") dXt = rXtdt + aXtdWt, Xo = 86.3, r = 0.03, a = 0.05. Abscissa axis: the optimal 
grids, Xtf. = xl, tk = fcA, A — 0.02, k — 1, . . . , 25, i = 1, . . . , Nk- Ordinate axis: the associated weights, f(Xt^. — x\), 
fc = 1, . . . , 25, i — 1, . . . , Nk- Xt^ is depicted in dots '•', Xt25 i^ represented by the symbol '*', ii — 0.02 and f25 ~ 0.5 and the 
remaining in continuous line 
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X,^, ^Xi, + t(X,,|a + it{X,J\/AZ, Z'-MGA); 6(3:) ^ ra-; iT{r) = i)- 




X,, = 4, ffc=itA,fc = l 25; A = 0,02;i = l, ,..,jV(. 



Figure 2: ("Pseudo-CEV model") dXt = rXtdi + i9(Xf+V(l + Xt2)-^''2)dWt, Xo = 100, r = 0.15, i? = 0.7, 5 = 0.5. 
Abscissa axis: the optimal grids, Xt^ = x]., tk = fcA, A — 0.02, fc = 1, . . . , 25, i = 1, . . . , A'^fc. Ordinate axis: the associated 
weights. Xti is depicted in dots '•', Xt25 '^ represented by the symbol '*', ii = 0.02 and f25 = 0.5 and the remaining in 
continuous line. 
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